RANDOM FIELD MODELS FOR RELAXOR FERROELECTRIC BEHAVIOR 
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Heat bath Monte Carlo simulations have been used to study a four-state clock model with a 
type of random field on simple cubic lattices. The model has the standard nonrandom two-spin 
exchange term with coupling energy J and a random field which consists of adding an energy D to 
one of the four spin states, chosen randomly at each site. This Ashkin- Teller-like model does not 
separate; the two random-field Ising model components are coupled. When D/ J — 3, the ground 
states of the model remain fully aligned. When D/J > 4, a different type of ground state is found, 
in which the occupation of two of the four spin states is close to 50%, and the other two are nearly 
absent. This means that one of the Ising components is almost completely ordered, while the other 
one has only short-range correlations. A large peak in the structure factor S{k) appears at small 
k for temperatures well above the transition to long-range order, and the appearance of this peak 
is associated with slow, "glassy" dynamics. The phase transition into the state where one Ising 
component is long-range ordered appears to be first order, but the latent heat is very small. 



I. INTRODUCTION 

Shortly after the seminal work of Imry and Ms} on the 
effects of random fields in ferromagnets, it was pointed 
out by Halperin and Varma^ that similar ideas could be 
used to understand the effects of atomic disorder in per- 
ovskite ferroelectrics. Halperin and Varma showed that 
atomic disorder which coupled linearly to the ferroelectric 
order parameter behaved in an equivalent fashion to the 
random field in a ferromagnet, and thus could cause the 
ferroelectric phase transition to become smeared. They 
seem to have been unaware, however, that ferroelectrics 
displaying such smeared transitions had been discovered 
long before 1^ These materials, now called relaxor fer- 
roelectrics, were originally referred to as ferroelectrics 
with diffuse phase transitions. The random- field ideas of 
Halperin and Varma did not immediately become popu- 
lar in the field. In the 1984 review of Isupov^ the diffuse 
phase transitions were modeled using only random-bond 
disorder. However, by the early 1990's it was generally 
recognized that random-bond disorder was not enough 
to explain all of the observed effects, and models in- 
cluding both random bonds and random fields became 
widely usedi^iSiL^ Until now it has been assumedZ*^ that 
in the presence of a strong cubic crystal-field anisotropy 
one could treat the different Ising components indepen- 
dently. In this work we will show that for a certain type 
of random field new phenomena occur which cannot be 
understood without including coupling between the com- 
ponents. 

The identification of relaxor ferroelectrics as a broad 
class of materials in which electric dipoles behave in 
a manner analogous in many respects to the magnetic 
dipoles in spin-glasses was made forcefully by Burns 
and DacoljS who generalized and extended ideas of 
CourtensiiS It is essential, however, to note that, as 
shown by Halperin and Varma^S in the case of electric 
dipoles it is almost inevitable that the disorder which is 



an essential feature of these materials will produce ran- 
dom fields. The presence of these random fields lim- 
its the development of long-range spin-glass order in 
three dimensionsiii Therefore, it is not expected that 
the glassy behavior which is seen in the disordered elec- 
tric dipole materials represents a true phase transition 
into a state with long-range spin-glass-type order. 

We would like to study the development of spatial cor- 
relations among the electric dipoles. In order to do this, 
we will work with a model which does not retain all of 
the atomic details. In the model studied here, we crudely 
represent each perovskite unit cell by only one dynami- 
cal variable. We will refer to this variable as a spin, al- 
though it actually represents atomic displacements in the 
unit cell. The model includes a cubic crystal anisotropy, 
a nearest-neighbor two-spin interaction, and a random 
field. 

As discussed by Fire and Blinc,^ if we assume that 
the cubic crystal anisotropy is so strong that each spin 
points along one of the [111] directions and the random 
field is represented by a simple vector coupling linearly 
to the spin, then the different components of the spin act 
independently.^^ Thus the model reduces to a set of de- 
coupled random-field Ising models. We may also include 
randomness in the bonds. (To perform this decoupling 
on Eq. (1), we use a coordinate system rotated by 45°.) 
This model has been studied in detail over a number of 
years. It is believed to be generally well-understood, al- 
though there remains some controversy about the values 
of the critical exponentsi-2*i^ 

Similar models have been studied in which the spin 
components do not act independently. Then the behav- 
ior does not reduce to that of an Ising model. For exam- 
ple, there are the three-state Fotts model with a random 
fieldji^ and the four-state Fotts model with a random 
fieldfi^ which are still qualitatively similar to the random- 
field Ising model. When the random field becomes strong 
enough, the ground state of the system breaks up into 



Imry-Ma domains, and the long-range order is destroyed. 

There is also the cubic model with random 
anisotropyr^ which, when the randomness is strong, 
shows a new type of domain phase with long-range or- 
der. This phase may be thought of as one in which each 
of the spin components orders independently, but differ- 
ent parts of the lattice are dominated by different com- 
ponents. Thus, in this phase, there is a network of 90° 
domain walls embedded in the ordered phase. A some- 
what similar domain state was found by Kartha, Castan, 
Krumhansl and Scthna^* in a model of shape-memory 
alloys. 

The existence of such a domain state with long-range 
order has been proposedi^iSli as the explanation of re- 
laxor ferroelectric behavior. Further, EgamiSi has em- 
phasized that the network of 90° domain walls is probably 
an essential element leading to the large piezoelectric re- 
sponse of relaxor ferroelectrics. Therefore, we would like 
to study a random-field model which has such a phase. 

II. THE MODEL 

The simplest model which has the desired properties is 
a four-state vector Potts model (four-state clock model) 
with a random Potts field, on a simple cubic lattice. The 
Hamiltonian of this model is 

H=-jY,^,-Sj+hrJ2Ss,-n,- (1) 

Each spin Sj is a dynamical variable which has four al- 
lowed states: (1,0), (0,1), (-1,0) and (0,-1). Each is an 
independent quenched random variable which assumes 
one of these four allowed states with equal probability. 
The (ij) indicates a sum over nearest neighbors of a sim- 
ple cubic lattice. The extension to three-component spins 
is conceptually straightforward, but will not be discussed 
in detail here. 

The sign of hr is chosen so that for large positive val- 
ues the probability of spin Si in state is strongly sup- 
pressed. Density-functional calculations for a typical re- 
laxor ferroelectric alloy^ have shown that the most sig- 
nificant interactions created by the alloy disorder are re- 
pulsive and short-ranged. Thus we believe that this is 
a reasonable first approximation to the interactions pro- 
duced by the disorder. 

When hr is negative, the random Potts field favors one 
state on each site, just as a vector field would. In that 
case there is no qualitative difference between the random 
Potts field and the random vector field. ^^-^"^ Thus we will 
not present any calculations for negative here. 

When hr is large and positive, we can approximate the 
random field term as a projection operator which forbids 
Si from occupying the state n^. Then within mean-field 
theory it becomes simple to calculate the approximate 
ground states. There are eight of these mean-field ground 
states. For example, in one such ground state 75% of the 
spins are in state (1,0), and the other 25% of the spins. 



which are prevented from being in this spin state by their 
local rii are in state (0,-1). The energy of this state for 
a simple cubic lattice is easily calculated to be —1.875 J 
per spin. Additional details of the mean-field theory are 
discussed in the Appendix. 

As we will see, the energies of the ground states found 
by computer calculation are less than —2 J, and they dif- 
fer qualitatively from the mean-field ground states. The 
problem with the mean-field ground states is easy to un- 
derstand. If we consider the above example, on a simple 
cubic lattice the 25% of spins in the state (0,-1) are bro- 
ken up into finite clusters. If we move all of the spins in 
such a finite cluster into the state (0,1), the total energy 
does not change. One might expect, therefore, that for 
large hr, the long-range order in the ground state would 
still be along one of the four [1,0] directions, since the 
25% of the spins in the finite clusters should not be able 
to exhibit long-range order. 

However, a true ground state on the simple cubic lat- 
tice actually has a very different structure. It turns out 
to be possible to find states with, for example, 50% of 
the spins in state (1,0) and 50% of the spins in state (0,- 
1) which have a much lower energy than the mean-field 
ground states. As we will show, to find states which look 
like the true ground states for large hr within mean-field 
theory, we need to add another term to Eq. (1). 

The author finds it helpful to compare the model stud- 
ied here to the results for the three dimensional Ashkin- 
Teller model,^^ whose Hamiltonian is 

Hat - J ^{o-i ■ o-j + • Tj) - 2 J4 ^((Ti • (Tj){Ti ■ r,) , 

(ij) {ij) 

(2) 

where <t = ±i(l,l) and t = ±i(l,— 1) are Ising vari- 
ables. (Note that our notation and normalization differ 
from those of Ditzian et alm^) 

Each of the three linearly independent components of a 
random field of a four-state model can now be identified 
with one of the terms of Eq. (2), and thus one of the 
mean-field order parameters. Expressing the random- 
field term using the a and r variables, and adding it to 
Eq. (2), we get 

Hrfat = i?AT+(/ir/4)^[l+2(cr,-n,)+2(Trn,)+477((Ti-n 

i 

(3) 

If we set 77 = 1 and J4 = 0, then Eq. (3) becomes just 
Eq. (1) expressed in the coordinate system defined by <t 
and T. If we then set rj = 0, the random Potts field term 
is replaced by the usual random vector field term. 

In a real alloy we would not expect that exactly one 
of the local states was blocked on each site. Instead, we 
would expect to find a random distribution, with some 
sites having no blocked states, some with one blocked 
state, etc. It would be straightforward to generalize the 
random field term to allow this, by letting 77 in Eq. (3) 
become a function of i. For example, we could use a 
probability distribution for r]i which had weight at ±1. 
Such a generalization is likely to be useful in modeling 
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the properties of specific alloy systems. However, a more 
general model has more parameters. Exploring its entire 
parameter space would be tedious, and seems premature 
at this time. 



III. PHASE DIAGRAM 

If hr/ J is chosen to be an integer, then the energy of 
any state is an integer multiple of J. Then it becomes 
straightforward to design a heat bath Monte Carlo com- 
puter program to study Eq. (1) which uses integer arith- 
metic for the energies, and a look-up table to calculate 
Boltzmann factors. This substantially improves the per- 
formance of the computer program, and was used for all 
the calculations reported here. (If desired, one could do 
almost as well for half-integer values of hr/J.) Lattices 
with periodic boundary conditions were used throughout. 
Two different linear congruential random number gener- 
ators were used: one to choose the random n,j and one to 
choose the Boltzmann factors for the Monte Carlo spin 
flips. 

A series oi LxLxL samples of various lattice sizes and 
values of hr was studied to map out the phase diagram, 
which is shown in Fig. 1. For hr/J < 3 the behavior of 
the random Potts field model remains generally similar 
to the behavior of the random vector field model. At 
low temperature the system develops long-range order, 
with the order parameter aligned along one of the spin- 
coordinate axes. It is not required, however, that the 
critical exponents for the phase transition here will be 
those of the random-field Ising model, as happens for the 
vector random-field case. 

The (<t) phase found by Ditzian, Banavar, Grest and 
KadanofP^ has many of the features of the domain FE 
phase which we find in our model for large hr and low T 
. The (cr) phase of Eq. 2 only occurs for J4 < 0. This 
indicates that one of the effects of the random Potts field 
in Eq. (1) is to generate a negative effective value of the 
J4 coupling. This effect was also seen in the random- 
anisotropy cubic model?!^ It is the random-anisotropy 
part (the 7] term of Eq. (3)) which generates the J4 ef- 
fective coupling under rescaling, and prevents the decou- 
pling of Eq. (1) into two Ising models. If we want to 
study a mean- field theory for Eq. (1) wc should include 
a J4 term, thus arriving at Eq. (3). The effective value 
of J4 will depend on hr- In this way we can obtain a 
mean-field theory which reproduces the phase diagram 
of Fig. 1. 

A necessary condition for the existence of stable long- 
range [1,1] domain FE order is that there be infinite con- 
nected clusters of both of the two dominant spin states. 
Since it becomes impossible to fulfill this condition when 
the amount of long-range order is small, the transition 
from the domain FE phase to the paraelectric phase must 
be first order. This, however, is specific to our model with 
only first neighbor interactions on the simple cubic lat- 
tice. It is possible that this phase transition could become 



continuous if, for example, second neighbor exchange is 
included in the Hamiltonian. 

As discussed by Ditzian et alm^, even in the absence of 
randomness the details of what happens as we move from 
the (cr) phase to the Baxter phase are not clear. There 
may be a small area in the region of the [1,0] FE to [1,1] 
domain FE phase boundary in which an intermediate 
phase is stable. This would be similar to the intermediate 
phase which has recently been found experimentallyi2& 
We have indicated the uncertainty about what is going 
on in this region of the phase diagram by the question 
mark in Fig. 1. A generalization of the (cr) phase also 
exists for the three-component Ashkin- Teller model in 
three dimensions i2i 

The dotted line in Fig. 1 represents the approximate 
location of the onset of thermal hysteresis. Above the 
dotted line, the system will quickly relax to a state which 
is independent of initial conditions, but below this line 
the system retains a memory of initial conditions for a 
long time. This line does not represent any true thermo- 
dynamic singularity. Its precise location depends some- 
what on the size of L and the length of the Monte Carlo 
runs. Thus the dotted line should be identified as the 
glass temperature, Tg. For small values of L the system 
can develop order on the scale of L above Tg. For L = 64 
and /ir > 3 it is not possible to find the ground states 
in a reasonable time by cooling the system from high 
temperature. 



IV. SPECIFIC HEAT 

In Fig. 2 we display the specific heat for Eq. (1) using 
four different values of hr- These curves were calculated 
by numerically differentiating and averaging Monte Carlo 
data for the energy of runs for four L — 6A lattices for 
each value of hr, starting at high temperature and cooling 
slowly. For each run, the temperature was changed in 
steps of 0.05J. 40,960 Monte Carlo passes through the 
lattice were performed at each value of T, with data being 
collected after each 20 passes. The first 248 data points 
at each temperature were discarded, and the remaining 
1,800 data points were averaged. A similar procedure 
was used for heating runs. 

Note that for hr — (not displayed), the specific heat 
for Eq. (1) reduces to a that of an ordinary Ising model on 
the simple cubic lattice, with Tc/ J — 2.256. This model 
has a very sharp singularity at Tc, whose height diverges 
to infinity with L. We see in Fig. 2 that as hr increases, 
the peak in the specific heat broadens and shifts to lower 
energy, with approximately half of the change occurring 
by hr = 3. By comparing with Fig. 1, it is seen that for 
hr > the peak in the specific heat is centered above 
the region of glassy behavior. The average heights and 
widths of the peaks shown in Fig. 2 are essentially the 
same for L — 32 (not shown) as for L = 64. As one would 
expect, there is more variation from sample to sample for 
the smaller lattices. 
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In Fig. 3 wc compare data from the cooling runs with 
hr = 6 with heating runs on the same lattices. The 
initial condition for each heating run was a state with 
long-range order. The temperature at which the long- 
range order collapses is identified as Tg. We see that 
below Tg there is a region of T where the specific heat is 
slightly higher for the cold-start sample. The crossing of 
the two curves near T/J= 1.2 gives us a crude estimate 
for the ordering temperature Tc. 

It is often possible to estimate Tc for a first-order phase 
transition with high precision. In order to do this with a 
computer simulation, however, it is necessary to run the 
simulation close to T^, with L larger than the correlation 
length, for a time which is long enough so that the sample 
can spontaneously order and disorder several times. We 
are far from being able to satisfy that condition here, due 
to the glassy dynamics for T <Tg. 

The integrated area between the cooling curve data 
and the heating curve data gives us an estimate of the 
value of the latent heat at T,. . Wc can also look directly at 
the difference in energy between the heating run and the 
cooling run at T = 1.2 J. Thus we find that the latent is 
about 0.01 J. The fact that the slope of the heating curve 
is so similar to that of the cooling curve and the latent 
heat at Tc is so small indicates that there is very little 
difference in the local structure of this model between the 
ordered states of the [1,1] domain ferroelectric phase and 
the metastable states which are found by slow cooling. 
A small latent heat is a natural consequence of a large 
(but finite) correlation length at Tc. 



V. CORRELATION FUNCTIONS 

In order to obtain information about the two-spin cor- 
relations, we calculate the angle-averaged structure fac- 
tor S{k). In Fig. 4 we display a log-log plot of S{k), 
obtained by averaging over the data from the cooling 
runs for the L = 64 lattices, as a function of T/J. We 
see that for all these values of /i,. > 3, the qualitative 
behavior of a sample cooled from high T without any ex- 
ternal ordering field is the same. As we lower T, the peak 
near k = first grows and then saturates. The rate of 
growth of this peak is maximal at the temperature where 
the specific heat has its maximum. As hr decreases, the 
height of the peak and the value of the correlation length 
increase. For hr — 3, the correlation length at low T 
becomes larger than our sample size. However, these 
zero-field-cooled samples never show true long-range or- 
der. With true long-range order, the size of the small-fc 
peak would decrease as T decreases below Tc, because 
the intensity in the long-range-order 5-function will not 
appear on the log- log plot. 

In Fig. 5 wc show data for S{k) from the heating nms, 
where the samples begin in an ordered state. Under these 
conditions, the size of the small-fc peak continues to grow 
as T is increased above Tc, until T reaches Tg and the 
sample can equilibrate. In Fig. 5(c) we show data from 



runs which are ten times longer than the standard runs. 
The hr = 4: samples arc able to equilibrate at T/J = 1.7 
during the long run {i.e. they become indistinguishable 
from the data from the corresponding cooling run), but 
the changes at T/J = 1.5 are small. By T/J = 1.4 
increasing the length of the run by a factor of ten has a 
negligible effect on the observed state of the sample. 



VI. GROUND STATES 

Somewhat surprisingly, it turned out to be possible to 
find approximate ground states of the hr = oo samples 
by a simulated annealing procedure. In Fig. 6 we show 
S{k) for the approximate ground states of the L — 64 
lattices with hr = oo, averaging over each of the four [1,1] 
directions for each of the four samples. The data shown 
in the figure arc fit quite well with by a Lorentzian line 
shape, with a correlation length of about 10 lattice units. 
This diffuse-scattering part of S{k) contains about 54% 
of the spectral weight. The other 46% of the weight is in 
the (5- function at A; = 0, which is due to the long-range 
[1,1] order. 

What this means in detail is that for large hr a ground 
state contains approximately 49% each of the two major- 
ity spin states, and about 1% each of the two minority 
states. The energy per spin of the lowest ground states 
found for each sample was about —2.010 J, with the en- 
ergies of the "groimd states" of the same sample in the 
other three [1,1] directions lying typically about 0.003 J 
higher. Since L = 64 is large compared to a correla- 
tion length of 10, this result for the ground state energy 
should represent the large L limit, while the differences in 
energies between alternate ground state directions of the 
same sample should scale like 1/ L^/"^. Note that these ac- 
tual ground state energies are substantially below the en- 
ergy of the simple mean-field ground state of the hr = oo 
model. 

For L = 64 and hr = oo the energy of a zero-field- 
coolcd state at T = 0, which has no long-range order, is 
only about 0.01 J higher than the energy of a true ground 
state. Thus for lattices smaller than about L = 32 the 
energy differences between the ground states of different 
[1,1] directions become larger than the energy difference 
between a zero-field-cooled state and a true ground state. 
Under these conditions the zero-field-cooled sample is 
able to reach a long-range-ordered state. The zero-field- 
cooled L = 64 samples can be decomposed into domains 
of different [1,1] ground states, with 90° walls between 
the domains. 

The probabilities of having effective fields of magnitude 
6, 5, ... and in a ground state, averaged over the four 
ground states for each of these L = 64: hr = oo samples 
are 0.1121, 0.2780, 0.3148, 0.1685, 0.0751, 0.0383, and 
0.0131, respectively. Thus, in a ground state only about 
one spin in nine is surrounded by fully aligned nearest 
neighbors. More than 1% of the spins are in zero effective 
field. 
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Some of the spins can flip freely between two positions 
in a ground state with no cost in energy. Since the spins 
in this model are discrete and two-dimensional, being 
in zero effective field is not a necessary condition for a 
spin to be "free" in this way. The bulk of these free 
spins are flipping between the two majority spin states 
for that particular ground state. It is easy to find that 
the fraction of free spins in a particular ground state is 
about 6.5%. However, a spin which is not free in one 
ground state may become free by the flipping of other 
free spins. The fraction of spins which can become free 
by the flipping of other free spins, one at a time, is 11.9%. 
Since 0.119 is much less than the critical concentration 
for uncorrelated site percolation on a simple cubic lattice, 
it is not surprising that a set of free spins deflned in this 
way consists of small isolated clusters. 

The probability that a free spin may occupy all three of 
its allowed spin states with no energy cost is negligible to 
this level of accuracy. Thus a lower bound on the ground 
state entropy per spin for large hr is 0.065 ln{2), and 
0.119 ln(2) is an upper bound. Calculating the ground 
state entropy of this model is a complex problem, because 
it requires flnding all of the ground states. 

Thus each of the four " [1,1] ground states" of a sample 
is actually a cluster of ground states in the phase space 
of the model. The simulated annealing procedure works 
because, if we use one of the simple mean-field ground 
states, which has 75% of the spins in one spin state and 
25% of the spins in a second spin state, as an initial con- 
dition, the direction of the order parameter will rapidly 
relax to the closest of the [1,1] energy minima, even at 
temperatures well below Tc- The fact that this is possible 
is partly due to the large number of spin rearrangements 
which can be made with no energy cost. This prevents 
the system from being easily trapped in a metastable 
minimum of the free energy which is close to the initial 
mean-field ground state. 

Although a ground state energy of —2.01 J means that 
about 33% of nearest neighbor spin pairs are oriented at 
90° from one another, the fraction of nearest neighbor 
spin pairs which are pointing in opposite directions in 
one of these ground states is only 0.00011. The small 
number of spins in the minority states exist in compact 
blobs whose diameter is approximately the correlation 
length of 10 lattice units. 

For hr/J = 3, every sample studied had four ground 
states, each one fully aligned along one of the four [1,0] 
spin states. There is an exponentially small probability 
of having an unusual local configuration of the which 
would misalign a small cluster of spins in a ground state, 
but this was never observed. Thus, ignoring the statisti- 
cal fluctuations in the sample average of n^, the ground 
state energy in the [1,0] ferroelectric phase is found to 
he Eq — ~ 3 J + 0.25 hr. The value of Eq for any 
value of hr cannot be greater than —2.01 J, its value for 
hr = oo. From this we expect that the boundary between 
the [1,0] phase and the [1,1] domain phase should be at 
hr/J slightly less than 4, as observed. 



VII. DISCUSSION 

In the past, various kinds of evidence have been pre- 
sented to argue that relaxor ferroelectrics represent a 
distinct class of materials, and are not merely the in- 
evitable consequence of adding some alloy disorder to 
any ferroelectric. De Yoreo, Pohl and Burns2& stud- 
ied the low-temperature properties of a variety of ferro- 
electric alloys, and found that they could be separated 
into two classes. Relaxor ferroelectrics showed glassy 
low temperature behavior and normal ferroelectrics did 
not, even when the normal ferroelectrics were random al- 
loys. Recently, Viehland et a^i ^?'?*^'?^ have argued that 
in some materials they find a well-defined transition be- 
tween a normal ferroelectric phase and a relaxor ferro- 
electric phase. It would be remarkable if the low tem- 
perature thermal properties of such a crystal could be 
switched from glassy behavior to crystalline behavior by 
poling the sample. To the author's knowledge, this ex- 
periment has not yet been tried. 

For the two-spin exchange interaction of Eq. (1), a 180° 
domain wall has twice the energy of a 90° domain wall. 
As a result of this, 180° domain walls become rare at low 
temperature in this model. As pointed out by Biirgel, 
Kleemann and Bianchi,— ferroelastic interactions raise 
the energy of 90° walls, but do not have much effect on 
the energy of 180° walls. 

To incorporate this effect into our model, we would add 
a biquadratic term to Eq. (1) of the form —K{Si ■ Sj)^. 
When K/J becomes large, the structure of the domain 
state will change^SS because 180° domain walls now have 
a lower energy than 90° domain walls. It is difficult to 
study this Biirgel, Kleemann, Bianchi domain state with 
computer simulations if the size of a ferroelastic domain 
is large. A sample which is the size of a single ferroelastic 
domain (or smaller) will behave essentially as a random- 
field Ising model. 

VIII. CONCLUSION 

In this work we have used Monte Carlo computer sim- 
ulations to study a simple model of a ferroelectric alloy. 
In this model, the alloy disorder causes some of the posi- 
tions of the polarizable ions to be strongly disfavored. We 
have found that this model displays a type of ferroelectric 
domain phase which does not exist in the model of Pirc 
and Blinci^ This phase may be thought of as resulting 
from the addition of a random field to the (cr) phase^^ of 
an Ashkin- Teller model. It seems likely that this model, 
and some natural generalizations of it, will help in the 
understanding of relaxor ferroelectric behavior. 

APPENDIX: MEAN FIELD THEORY 

The presence of the random field greatly complicates 
the mean- field theory for Eq. (1). It is necessary to distin- 
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guish the four classes of spins corresponding to the four 
different allowed values of n. Label the four directions 
(1,0), (0,1), (-1,0) and (0,-1) as direction 1, 2, 3 and 4, 
respectively. Define to be the probability that a spin 
of class ^ points in the v direction. Then ^^p^i^ — 1 
for any fi. But, in the general case, we have to deal with 
twelve independent variables in the mean-field theory. If 
we take the limit that hr becomes infinite, there are only 
eight independent variables, because we can then set p^^^, 
to when ^ = v. 

Let z be the number of neighbors of any spin. Then 
the mean-field expression for the energy per spin is 

S=-(zJ/32) Vi..Vp,'v'C0s{T:{v~v')l2). (A.l) 

In addition to the usual factor of 1/2 to avoid double- 
counting, there is factor of 1/16 arising from the assump- 
tion that each value of \x (i.e. each n) is equally weighted. 



The entropy per spin is 

5 = -l/4^p^,ln(p^,), (A.2) 

and the free energy is, as usual, F — E — TS. Then solv- 
ing the mean-field theory requires finding the minimum 
of F{T) in the phase space of the Pf^^. 
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FIG. 1: Phase diagram of the random Potts field model on 
simple cubic lattices, showing the paraelectric (PE), ferro- 
electric (FE) , and domain EE phases. The solid lines indicate 
first-order transitions, and the dashed line indicate the ap- 
proximate onset of glassy dynamics in the PE phase. The 
question mark is discussed in the text. 



EIG. 3: Specific heat vs. temperature for the random Potts 
field model on 64 x 64 x 64 simple cubic lattices, for hr = 
6, averaging data from four samples. The solid line shows 

zero-ficld-cooled data using a random initial condition, and 
the dashed line shows zero-field-heated data using an ordered 
initial condition. 



FIG. 2: Specific heat vs. temperature for the random Potts 
field model on 64 x 64 x 64 simple cubic lattices, for various 
values of hr- Data from zero- field-cooled runs, averaging data 
from four samples. 
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FIG. 4: Angle-averaged magnetic structure factor at a se- 
quence of temperatures for the random Potts field model on 
64 X 64 X 64 simple cubic lattices, log-log plot. The points 
show averaged data from four samples, using zero-field-cooling 
and a random initial condition, (a) hr — oo; (b) hr = 6; (c) 
hr = 4; (d) hr = 3. 



FIG. 5: Angle-averaged magnetic structure factor at a se- 
quence of temperatures for the random Potts field model on 
64x64x64 simple cubic lattices, log- log plot. The points show 
averaged data from four samples, using zero-field-heating and 
an ordered initial condition, (a) hr — 6; (b) hr = 4; (c) 
hr — 4, relaxed (see text); (d) hr = 3. 



FIG. 6: Angle- averaged magnetic structure factor at T = 
for the random Potts field model with /ir = cxd on 64 x 64 x 64 
simple cubic lattices, log-log plot. The points show averaged 
data from four samples, using four states from each sample, 
one from each [1,1] direction. 
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